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ABSTRACT 


This thesis develops a means by which a closed loop frequency 
response may be calculated by digital computer, given only the circuit 
configuration and component transfer functions. Using these data, a 
program is developed to perform the necessary logic for circuit analysis 
and obtain a performance function and frequency response over any de- 
sired frequency range, for an input and output at any specified circuit 
location. 


The method is demonstrated by obtaining the performance func- 
tion and frequency response of a complicated closed loop system. 
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OBJECT 


The object of this work is to determine the feasibility and 
devise a method for numerical computation of closed loop fre- 
quency response and to apply the method to an illustrative pro- 


blem. 
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CHAPTER 1 


INTRODUCTION 


In the design and analysis of feedback control systems it is 
often necessary or deSirable to know the steady state sinusoidal 
frequency response of the system. Various methods for obtaining 
this frequency response are available to the engineer, the most 
often used, perhaps, being those based on the Bode Diagram and 
Nichols Chart. For complicated high order systems, these graph- 
ical techniques are apt to become tedious and time-consuming; even 
more disadvantageous is that the time must be spent by the engineer, 
or at best, an assistant with knowledge of the engineering principles 
involved. There is also to be considered the human factor and its 
accompanying probability of error. | 

It would be advantageous then to have available a method re- 
quiring a minimum application by the engineer, one that is relatively 
fast, free from the human factor and, after the minimum initial ap- 
plication, could be handled by less-technically-qualified personnel, 
With the advent of high speed, large data capacity computers, it is 
reasonable to assume that such a method could be made available. 

It is toward this end that the method to be presented in this paper is 
directed. 

In general, the minimum input information required from the 
design engineer is considered to be: 

1. the transfer functions of the individual components of 
the system 
a means of describing how the components are connected 
the frequency range and desired increment of frequency 
over this range for which the response is required. 
Item 1 is assumed known to the engineer and need only be entered in 


the form described in Chapter 2. Item 2 is desired to be as simple as 


ae 





possible but still be of such a nature that a poSitive check may be made 
by the engineer to enSure that the information is entered correctly. 
This is described in detail in Chapter 2. Item 3 consists only ofa 
statement of the values noted. 

It is noted here that no new mathematical concepts are Set forth, 
in fact, the simplest of arithmetic operations are involved. These basic 
fundamentals are the most readily adapted to digital computation. 
Graphical solutions such as the Bode Diagram were designed not be- 
cause the fundamental mathematics were unknown, but because these 
Simple arithmetic operations were so laborious as to be prohibitive. 
The modern computer performs these operations at an extremely 
rapid rate and permits programming the heretofore prohibitive methods. 

Chapter 2 contains the physical interpretation, with appropriate 
block diagrams and flow charts, of how the program works; i.e., how 
the information is read into the computer and the general steps neces- 
sary to compute the frequency response. The physical reasoning and 
diagrams are intended to be general in that they describe the necessary 
programming steps and hence, are applicable to any digital computer. 
At the same time, they are in great enough detail to permit the pro- 
gram to be adapted to any computer with a minimum amount of additional 
labor. 

Chapter 3 and appendices thereto contains the detailed program 
as written for the IBM 650 Computer and its associated equipment RAMAC, 
The "MAC" language in use at the Instrumentation Laboratory at MIT is 
utilized because of the relative ease of programming. 

Chapter 4 contains a brief summary of the salient features of the 
program and its application to some sample systems, 

The primary objective of this thesis is to set forth a general pro- 
gram, stressing Simplicity of use and not necessarily to derive the opti- 
mum program with respect to minimum time and compactness, For this 
reason, it is realized that the program of Chapter 3 may not be the ulti- 
mate in computer programming but is general enough to be adapted to 


any modern digital computer. 





Crear mene 2 


GENERAL DESCRIPTION OF THE SYSTEM 


pel. General 

The program that is used to determine a general form of the 
system performance function and from this compute the frequency res- 
ponse of a linear, time invariant system is described in this chapter. 
The description is meant to be quite general in that is applies to pro- 
gramming any digital computer. If not general, the program would 
be of little value for any computer except the specific one used here. 
Even though general, a certain amount of detail is necessary to en- 
able a prospective user to program the system without an unduly great 
amount of effort. The description in this chapter is taken by steps, 
beginning with a broad over-all aspect and block diagram flow chart 
and working down to the details. This gives the reader an opportunity 
to decide just how far he wants to go and permits comprehension at 


any desired level of detail. 


2.2. Overall Description 

As stated in Chapter 1, there are numerous methods in the litera- 
ture and in current usage to calculate the frequency response. Some 
of these, particularly the graphical techniques, are not readily adapt- 
able to digital computation. We may note, however, that in obtaining 
the system performance function by any method, we get at some Stage 


of the calculation an equation of the form: 
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where the G's represent the transfer functions of the system components. 
A more careful study reveals that the components of the denominator 
terms are thoSe components making up a loop or combinations of loops 
in the physical system. Although not quite so obvious, it turns out that 
some of the terms of the numerator may be recognized as loops and some 
as direct paths from the input to the output of the system. Therefore, it 
would seem that if we could determine the general scheme of how these 
loops and paths combine in the equation, then it might be possible to 
readily adapt this to digital computation, Mason (Ref. 5 ) obtains the 
general format of how they combine and arrives at a general equation 
representing the performance function for an input at i and an output 

at j as: 
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where: 
P is the product of transfer functions of all components in 
a direct path fromitoj. This product is defined as the 
"Path function" of path 'k"’. 
LE is the product of transfer functions of all components of 
a loop and is defined as the "loop function" of loop "r"'. 
L._ is products taken 2, 3, ... at atime of loop functions 
of all loops not "touching" one another. 


Ly, is the loop function of a loop not "touching" path "k" 


is the products taken 2, 3, ... at atime of the loop 
functions of all loops neither "touching" one another nor 
touching path '"'k" 

The only rule involved is that of signs. This may be stated as: "In trav- 
ers\ing a loop, if the number of (-) signs encountered is odd, the sign of 
the term is (+). If the number is even (including zero), then the sign of 
the term is (-)", This is readily understood from basic feedback system 


theory where we are familiar with the result that a negative feedback 





function leads to a characteristic eqvation 1 + the loop function. In 
traversing a ''direct path" if the number of (-) signs encountered is 
even (including zero) the sign of the term is (+); if odd, the sign of 
the term is (-). 

To demonstrate the use of the basic equation and the "touching 
vs. non-touching" character of loops and paths consider the following 


system where Gris the transfer function of block (N) 


Fin ¥ 






From the block diagram, it is seen that there are: 
1) 4 loops with loop functions L, = G,G,G,G,G 
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2 L,. = G,G,G,G4G_ + GeG,G_ + GoGr - G2G,Ge 
2) 2 direct paths with path 
functions P, = G, GoG3Gy4 
Dee NS 





The ''touching" of 2 loops can also be seen by inspection but if in doubt 
the question might be asked ''can a block be inserted in one of the loops 
without changing the loop function of the other?" If so, then they do not 


touch, In the block diagram above: 


loop (1) touches the remaining three 
loop (2) does not touch loop (3) but does touch (1) and (4) 


loop (3) and (4) do not touch one another 


then: 
Path (1) is seen to touch all loops: yD Loi = 0 
Path (2) touches all except loop (3): Z Lg = GoGy 


The performance function then is: 


+.G,. Gy. GaG 
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Ua G, Go G,G,G_.tG,G,GetGoGr ~ G,G,G_gtG, G,GgGoGr = GoGrG,G,Ge 


To obtain the frequency response, two courses are now available. 
The cbmponent transfer functions might now be entered and a program 


written to expand the performance function above into the form 


Asti pst™iicgt ee... 
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Then let S = j w and compute the desired response. This method would 
require extensive working and storage space and might well exceed the 


capacity of some computers. The other course available is to enter 





the transfer functions of the components, let S = j} w and compute a 
magnitude ratio and phase angle for each component for a particular 
value of w. Then, substitute these values into the au equation, 
perform the necessary mathematical operations and obtain the res- 
ponse for the particular value of w&. The program described in 
Chapter 3, uses the latter method. 

The form in which the final answer is received will vary with 
the computer; perhaps punched cards, or if a plotter is available, 
then the final result may be read directly from a plot. A flow chart 


of this broad description is contained in Fig. 2.1. 


2.3. Input Data 

It is desired that this be in as simple a form as possible and 
provide some sort of a check to the programmer that the informa- 
tion (data) has been correctly entered. To fulfill both of these re- 
quirements, it is considered best to use an arbitrary numbering 
system of all components and juction points (summers). Then the 
inputs to a numbered unit are the numbers of the unit (s) from which 
they came. All that is necessary to be entered then is the number of 
the unit and the number (s) of the input. To permit operation by the 
computer, these data are entered in an array. For example, con- 


sider the simple system illustrated here in block diagram form: 





If computer space is critical or if the system is extremely large, the 
blocks 2 and 3 may be just as well combined into a single unit. If con- 


sidered as two separate units, the input array would be: 


Unit Inputs 
2 4 
3 2 
a 1 -3 (the 1 is unnecessary but used here to 


illustrate the form for multiple inputs), 
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Note the (-) indicating this is to be subtracted in the summer; i.e. neg- 
ative feedback. 

The only other input data required are: input and output location num- 
bers, in this case 4,3; number of components and summers and mak- 
imum number of inputs to any unit. The reason for the first require- 
ment iS obvious; the last three are used to reduce computation time 
and conserve computer storage space by establishing the exact size of 
the array. 

For reasons to be explained later, it is necessary to number the 
summers also and is convenient to assign larger numbers to these 
than the highest numbered component and list them last in the array. 
The particular order in which either components or Summers are 


numbered is otherwise arbitrary. 


2.4. Computation of Performance Function 

Using the method described in paragraph 2.2., it is seen that 
what is actually required is a method of tracing through the numerical 
input data to determine the existent loops and direct paths of the sys- 
tem. These are found and recorded as groups of numbers; hence to 
determine the ''non-touching" character, it is necessary only to check 
numbers of each group against the other groups. Since it is possible 
for paths or loops to ''touch'' without having common components, the 
necessity for numbering the Summers is now Seen. It is not possible 
for loops to touch without having at least one common junction point 
(Summer) and vice versa. The broad scheme of this operation is 
described in the flow chart of Fig. 2.2. 

After determining the existent loops, paths, non-touching loops, 
etc., the groups of numbers are compiled to form the general expres- 
sion for the performance function, Since the requirement might exist 
for knowing the performance function, a Subroutine is included in the 
main program for compiling and punching out this function, before 
commencing calculation of frequency response. 

The methods used to determine the paths and loops are perhaps 
the heart of the entire program. The two methods are similar in 


general hence only the loop routine will be explained in detail and the 
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path routine in general except where it differs. The first step of the 
loop routine is to choose a Starting point in the array. Since any loop 
of the system must contain a junction point (summer), then we choose 
as the Starting point the first summer listed in the input data array 
(not necessarily the first one in the physical system). We could start 
with a component but can save time by considering only the summers. 
Hence, it is convenient to list the components first then the Summers, 
This order might well be reversed but the program of Chapter 3 is 
written for the order noted. The first summer is located by knowing 
the number of components other that summers. To enable computa- 
tion of the frequency response or in punching the performance function, 
the Summer numbers are not desired, This is obviated by assigning 
higher numbers to the summers than to the components, and discard- 
ing all numbers greater than that of the highest numbered component. 
The computer, after noting the starting number, then Scans the input 
array by rows until this number is found as an input to another num- 
bered unit. Then the new number is used as the "search number" 
until it is found as the input to another numbered unit, and so on until 
we have a complete loop. It is seen now that difficulty might be en- 
countered due to the output of a unit going to more than one other unit. 
This must be considered inthe program. The way this is done is to 
record the first output then continue the search to see if there are any 
more. If more are found they too are recorded at this time. It is 
impractical to find one complete loop then go back to the starting point 
and search for more loops. Hence, all points in all branches are re- 
corded at this time. This is more clearly described by the flow chart 
of Fig. 2.3. After recording (in the computer) all the loops contain- 
ing the first summer, the program continues to the second and re- 
peats the search and record process, ensuring that loops already re- 
corded are not repeated, and so on until all loops are found, It is 
realized, of course, that checks must be made for such things as 
traversing a loop only once and recording only complete loops. This 


is also brought out more clearly in Fig. 2.3. 
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Fig. 2.3. Detailed Flow Chart of Loop Routine 
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The path routine differs from the loop routine primarily in that 
the starting point is the number inserted in the input data as the sys- 
tem input for which the response is desired, and also the end of the 
path is the listed output number. The similarities and differences 
are seen more clearly by comparing the flow chart for the path rou- 
tine (Fig. 2.4.) with Fig. 2.3. 

The importance of these two subroutines and the program 
sophistication necessary to accomplish the desired result cannot 
be emphasized too strongly. The flow charts of Figs. 2.3. and 2.4. 
clarify the procedure satisfactorily for understanding the general 
method employed. However, they are not considered to be in great 
enough detail to permit writing a similar program for another com- 
puter. For this reason, the method of this subroutine is explained 
in great detail in Chapter 3. 

Having recorded all the existent loops and direct paths of the 
system, it remains to determine the remaining factors of the general 
equation, The idea is very simple; we merely check numbers of 
groups already recorded against appropriate other groups. If no 
identical numbers are found then the loops, etc., do not touch. Since 
this procedure is fairly readily understood, a general flow chart 
(Fig. 2.5.) is considered adequate. 

To obtain the performance function, we now combine the ap- 
propriate previously stored numbers and punch out the function term 
by term. Fig. 2.6., is a general flow chart of this routine. It 
should be pointed out that since it is foreseeable that the performance 
function might not actually be desired, then writing this subroutine 
is optional. For this reason a "by-pass" is put into the program and 
effected by entering a value of zero for a variable read in with the 
input data. If the performance function is desired, then any number 


other than zero is entered. 


2.5. Component Transfer Functions and Frequency Response 

Any component whoSe input-output performance can be deScribed 
by a linear ordinary differential equation with real constant coefficients 
has a transfer function which is rational in the Laplace variable S 


with the numerator and denominator polynomials expressible as 
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products of terms of the form A s? 7.0 + Can whichvany sof the 
coefficients may be zero. These transfer functions must now be 
entered. The form in which they are entered is discussed in detail 
with the program of Chapter 3 and illustrated in Chapter 4. 

Next a magnitude ratio.and phase angle is computed for this 
component by letting S = j] w and substituting the particular value 
of w for which we are computing the response. The required bits 
of information regarding frequency are: frequency range (upper 
and lower limit) and increment. (See Fig. 2.7.). 

To facilitate ease of recording the required information, a 
form is included with the program details of Chapter 3 and the 


sample system of Chapter 4, 


2.6. System Magnitude Ratio and Phase Angle 

Having computed, for the particular value of w, the magni- 
tude and phase angle of each component the next and final step is to 
substitute these values into the general performance function. If 
the entire performance function has been retained in storage, if 
and when punched, it could well be used to compute the system res- 
ponse; but, since it is desired to have this performance function 
calculation optional, then the routine for finding system response is 
written uSing the loops, paths etc., previously computed. This is a 
Straight-forward routine easily written and understood once it is 
known where the information is stored from previous Subroutines. 
Since the general equation from paragraph 2.2., contains only loops 
and paths and combinations of these, it is expedient to compute and 
Store first the magnitude and phase for each loop and path, then use 
these stored values in the combinations. The performance function 
is in general a summation of loops, paths, etc., hence a conversion 
from polar to rectangular form is also required and, after compil- 
ing, a conversion back to polar form. For any greater detail, the 
reader is referred to the program detail of Chapter 3 and the sample 


problem of Chapter 4. 
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2.7. Frequency Response "Re-Run" 

After obtaining a frequency response over an indicated range 
with a specified increment, it might be desired to change the range, 
increment or even the transfer function of some of the existent com- 
ponents. In this case, it would not be desired to recompute the gen- 
eral performance equation of paragraph 2.2., hence, provisions 
are made to re-enter the program with these different values. The 
method will vary with the computer used and for this reason, is not 
detailed here. A detailed step by step procedure for the IBM 650 


is contained in paragraph 3.10. 
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CHAPTER 3 


DETAILED DESCRIPTION OF PROGRAM WRITTEN FOR IBM-650 


ot. General 

It was stated in Chapter 1, that the MAC language (Ref. 6) is 
used in this program for simplicity of construction. Very briefly, 
MAC may be described as follows: A plug board wiring scheme per- 
mits recognition of certain special plain language commands such as 
"PUNCH", "GO TO", "IF", "DO" and "DO TO", Arithmetic symbols 
are also recognized; examples of these are the signs for "equals", 
"plus", 'minus", ''decimal point" and ''parenthesis'’. Any number of 
parentheses may be used consistent with proper arithmetical usage. 

Variables are created by several means, the most common one 
being an equation with the variable to the left of the equal Sign, e.g., 
A=2.5-+ B/2 is interpreted as; ''a new variable (or new value of a 
variable in use) is equal to 2.5 plus (B (a previously defined quantity) 
divided by 2). Upon encountering this command in a program, the 
value of A is computed and stored in the memory and may be uSed at 
any time in the right hand side of an equation. 

Routines for certain other functions, Such as trigonometric 
functions, Square root, absolute value, integer are built into the basic 
Storage and are referred to upon recognition of the functions: SIN, 
COS, TAN, ARCSIN, SQRT, ABS, INTEGER and others for which the 
reader is referred to Ref. 6 of the bibliography. Any number of sub- 
Scripts, either numbers or letters for which a value has been compu- 
ted, may also be used. It is noted that subscripted variables are treat- 
ed in the Same manner as non-subscripted ones and new variables may 
be introduced accordingly. If the subscript is itself a variable, its 


value muSt be defined. The subscript is then an "index" and muSt be 
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so indicated on the first card of the program, Variables are limited 
to letters of the alphabet, numbers, and combinations of either or 
both. The physical limits imposed by MAC is a total of 800 variables. 

The individual lines and cards of the MAC program are number- 
ed in Sequence, using the prefixes E, M, S to indicate exponent, main, 
and subscript cards respectively, Equation numbers when required 
appear to the right of the MAC statements and equations, and are used 
as a means of transferring control to selected points in the program 
or to establish loops. 

After all MAC program cards are read into the computer, an 
internal assembly routine converts this program into the basic language 
of the computer, begins at a designated point, and proceeds with the 
computation. This 'new program" may be punched on cards as it is 
being written, the cards then being referred to as a "load deck"; or 
if no future need is anticipated, the load deck is not punched. Use of 
this load deck in lieu of MAC cards for future operations reduces com- 


putation time considerably. 


3.2. Major Subroutines of the Program 
Reference to the basic equation for determining the performance 
function, Paragraph 2.2., shows that the component numbers are re- 
quired for: 
a. All closed loops of the system (L) 
b. All direct paths from input to output (P,) 
c, All combinations of 2, 3, 4 etc. loops not touching one 
another (L.. products) 
d. All loops not touching each path (Ly) 
All combinations of 2, 3, 4 etc. loops not touching each 
other and not touching each path (Loy products) 
In addition to the above, other major subroutines are: 
f. Calculation and punching of the performance function 
g. Calculation of frequency response (for each w) of: 
1. All components 
2. All loop functions, all path functions, certain pro- 
ducts of loop and path functions 
Denominator of the system response 
4. Numerator of the system response 


The system 
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Flow charts for each of these have been included in the general dis- 
cussion of Chapter 2. The complete program is appended as Appen- 
dix B and is discussed in detail in the following paragraphs. Each © 
major Subroutine will be explained separately, with reference to the 
appropriate flow chart and sections of the program. The numbers in 
parentheses following the subroutine paragraph title indicate the appro- 


priate location in the program. 


3.3. Input Data 
Before discussing the program, it is essential that the "Input 
Data Array" be understood. This is most readily accomplished by 


means of a block diagram of a hypothetical system: 


output 





Fig. 3.1. 


The numbering of components and summers is shown to be purely 
arbitrary except that summers are assigned higher numbers. Data 


Cards for MAC contain six or less numbers and are written in floating 
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point form (Ref. 7). The first card contains in the order listed: 


Maximum no, of inputs to any unit 


(excluding the system input) 2 
Number of component boxes 2 
Number of summers 3 
Input unit number if 
Output unit number 2 
Is performance function desired O for No. 

1 for Yes 


On the succeeding cards, the inputs to each unit are listed, using a 
separate card for each unit. Again the order of listing is arbitrary 


except that the summers should follow the component boxes: 


Unit No. Inputs to this unit 

1 8 0 0 (0-3) (numbers in parentheses 
2 6 0 0 (4-7) are not entered, but are 
3 7 0 0 (8-11) used for explanation fol- 
4 7 0 0 (12-15) lowing) 

3 2 0 0 (16-19) 

6 1 5 0 (20-23) 

7 -2 -1 0 (24-27) 

8 -5 0 (28-31) 

0 0 0 (3i2- 3a) 


The zeros are required since with MAC a "Read" command is 
used to enter data and this Read command must specify what values and 
how many are to be read. The program for which this Chapter is writ- 
ten uses the ‘maximum number of inputs to any unit" as the basis for 
specifying the width of the array and hence, the number of values to be 
read. The zeros are used as "stops" in the search routine. Note, the 
(-) signs indicating negative feedback into unit 7 and 8. 

For any system, a rectangular array is obtained, the width of 
which is the ''maximum number of inputs to any unit" + 2. The length 
of the array is "the number of component boxes + the number of sum- 


mers + 1 (the last row is all zeros). 
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Each position in the array is identified as a subscripted variable 
and the values read from the data cards are thereafter recognized as 
the value of the variable. The numbers in parenthesis are the sub- 
Scripts identifying the locations. Gis chosen as the variable to be 


subscripted. For example, the value of Goo Tome 


3.4, L.. Subroutine (M 0075 - M 0585) 

It is desired to determine the units comprising each of the closed 
loops of the system. As noted in Chapter 2, this is a very important 
but complicated subroutine. To ensure complete understanding requires 
that the program be explained in great detail. Since this is, of neces- 
Sity, a lengthy discussion the program details are included as Appendix 
A and only the method will be discussed here. Refer to the flow chart 
of Fig. 2.3. The basic method is to chooSe as a Starting point in the 
array, the first summer listed. Using this number as a "search num- 
ber” the input array is scanned, starting with G, until the number is 
found as an input to another unit. The Search number and the number 
of the output uhit are entered temporarily in an array designated by the 
Subscripted variable XX, The search continues for more outputs of 
the "Search number". If any more are found, a new line is formed in 
the XX~-array for each. After all outputs are recorded, the last output 
iS chosen as a new search number and the same procedure carried out 
again, writing a new line in the XX-array each time a new unit is found 
to be a point in the loop. Although not too obvious at this time, this is 
necessary to enSure the same path is not being traversed around a loop, 
and to have all numbers available when the loop is found to have been 
completed. 

To illustrate briefly, the first summer listed in the array of 


paragraph 3.3., is (6) which goes to (2), so the XX-array is: 


line (a) 6 : 
line (b) 6 2 


No more outputs of 6 are found so using (2) as a search number, its 


output is found to go to (5) and negatively to (7). These are recorded 
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in the XX array as lines (c) and (d): 


line (c) be “a2. 
line (d) -6 & 7 


(7) is the next Search number, for which an output is found to 3 and 4. 
These are recorded as lines (e) and (f). Then (5) is selected as a 


search number, the output of which feeds negatively to (8), (line g) 


limemte) 26 2 7 3 
linew(f) 208.2 7 4 
line (g) -6 2 5 8 


The routine continues, checking each new number to See if it forms 
a completed loop. If so, the appropriate line of the XX-array is re- 
corded permanently as a loop, for which the subscripted variable ZZ 
is used. Continuing the search process yields one new line from (g) 
and one from (f) but when (3) from line (e) is checked, it is found to 
go to (6) which is the starting point, hence a complete loop has been 
found and a new line is not entered in the XX-array, but is recorded 
as LL, to ZZ, with values of -6, 2, 7, 3 respectively. ZZ, is then 
recorded as zero aS a Separator from the next loop. 

After each number is located, it must be checked to See that it 
has not already been entered as a point in the line of the XX-array. 
If So, this represents a minor loop within the main loop that is being 
traced, This minor loop will be recorded as a line in a succeeding 


array. Hence, no new line is entered. This is illustrated as follows: 





line (h) “6928 8S @ 1 
line (i) —By "2? a " Se 
line (j) -G 2a. 4 88. 1 


Now use (1) from line (h) as a search number and find it goes to (6), 


a complete loop, and back negatively to (7). A new loop ZZ, to ZL,9 


Oe 





is recorded as -6, 2, 5, 8, 1 and a new line (k) is recorded as: 


line (k) +6 2 5 8 1 7 (note the sign change) 
_ line (1) Shue 1 ae 
line (m) 6 2 1 


Lines (1) and (m) are recorded from line (k) and we move back to the (1) 
of line (j). (1) is found to go to (6), a complete loop, hence, is record- 
ed in the ZZ list but not in the XX-array, and to (7); but since (7) has 
already been recorded in this line, it is disregarded and nothing is re- 
corded in the XX-array. The process continues until all loops of which 
summer (6) is a part, have been recorded in the ZZ list. The next 
summer of the "input data array'' is then chosen as a starting point and 
the Same procedure is carried out as for the first Summer, It is realiz- 
ed now that loops that have already been recorded from summer (6) must 
not be recorded again, so appropriate checks must be made to prevent 
this. When all summers have been checked for new loops, all of the 
existent loops of the system have been recorded and the program moves 
on to the next major subroutine. For greater detail as to how the above 


methods are programmed, refer to Appendix A. 


3.5, Subroutine for Determining Direct Paths (PK routine; MO590-1040) 

This routine differs only slightly from the LR routine above, and 
will be discussed in detail only where major differences occur. The 
search of the G-array is the Same as above except the idea now is to find 
all components in a path between two specified units. Therefore, the 
first number of the XX-array is the designated input unit and the array 
is written only once. 

Recording new lines in the XX-array (MO790-0820) and the checks 
made prior to this (MO735-0785) are the same as before except'here, the 
check for a complete path is made against the specified output and the 
duplicate point check includes the starting number. When complete paths 
are found, all units of the path are recorded as the subscripted variable 
YY and this procedure differs from the LR routine only in that the last 


number found is also recorded (MO970), 
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It is necessary to know the number of paths just as it was neces- 
sary to know the number of loops, hence a tally is made each time a 
new path is recorded (M1020 and M1040). Note that Yj is the number 
of paths and the subscripted variable Y also is used in the XX-~array. 
The reasoning here is that the variables used in writing the XX-array 
are not needed after this Subroutine is completed and with limited com- 
puter Space it is best to use variables that may also be used later. 
Throughout the program, the single letter variables are used many 
times but space is reserved for a designated number of variables with 
indexed subscripts (i.e., subscripts that are also variables) and, once 
used, this space is wasted unless the subscripted variables are chosen 
such that they may be used again. Since, as will be seen later, there 
are a number of constants associated with each path, (e.g. the number 
of loops not touching each path), it is convenient to use an indexed- 
subScripted variable to denote these constants. The subScripted vari- 
able Y is used for this and hence that space used by the Y's in writing 
the XX-array is not wasted. The only possibility of an overlap is Yo: 
Therefore, the tally of paths is kept as a single letter variable (M1020) 


and upon completion Yo is set equal to this tally. 


3.6. Products of Loop Functions Subroutine (M1050-1795) 

From the basic performance function equation of Chapter 2, it is 
seen that both the numerator and denominator are made up of a Series 
of terms, the number of which is determined by the complexity of the 
physical system. It is impractical to program an unlimited number of 
product terms and the question quite naturally arises as to how many 
Should the program be capable of handling. Intuitively, it is seen that 
the degree of complexity of the physical system also determines the 
amount of computer space required for the "input data array", the 
XX-array, the ZZ and YY lists, etc. There is therefore, a physical 
limit imposed on any program written in MAC without resorting to 
other means for obtaining computer space. Consideration of these 
limiting factors leads to an estimation of the maximum number of pro- 
duct: terms that should be programmed, The program appended to this 
report is written for a maximum of "products of loop functions taken 


four at a time of all loops not touching one another", 
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The method used to determine the touching vs. non-touching 
character of the loops is to check each number (component or summer) 
recorded in a loop against each number of every other loop. If no 
identical numbers are found in the two loops, then they do not "touch" 
and are to be So recorded. The details of the program are most 
readily explained by means of an example. 

Consider a system containing four loops. These have been re- 

corded in the ZZ list with a zero Separating the numbers of one 


loop from the next: 
) (2) (3) (4) 


LZ, LZ, LZ, LZ, 


1 at+2 b+ 2 C1 2 
a a) : : 
Za, AZ, Did, LD, 
ZZ, 4470 42499 Slog 2 ee ee 


Check loop (1) against loop (2) by first checking ZZ, against ZZ, 
to ZZ 


+ 2 


Bt 1° 
For convenience of locating loops by number at a later time, 


the first ZZ number of each loop is recorded as the check proceeds. 
The subscripted variable G is used for this since at the present time, 
the Space reserved for the G's is unused. Since a + 2 is not yet known, 
it is found by checking the ZZ list until a zero is found (M1135-1140). 
The next ZZ number after the zero is ZZ. +2 and Go is recorded as 
this value (M1160). Now the check consists of: ZZ, (first check 
point) against all numbers of loop 2. When ZZ, a9 (= 0) is found, 
pick up ZZ, as the "check point" and start again with ZZ. 4.5 a8 the 
"check number", The procedure is continued until 2D. , i) as 
found as the check point or until the check number and check point are 
identical. Since algebraic Signs are of no Significance, absolute values 
are used in the check. | 

If the check point becomes zero, loops (1) and (2) are recorded 


as a pair of non-touching loops, and loop (1) is then checked against 
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loop (3) and so on. When the loop number of the loop to be checked 
exceeds Zo (number of loops) (M1150) the procedure is to check loop 
(2) against the remaining loops in the same manner as before, except 
now the first ZZ number of each loop is known and uSed as the first 
"check number". The procedure continues until loop number (Zy - 1) 
is checked against loop number (Z,). 

The non-touching loop pairs found in the above process are 
recorded as the subscripted variable XX, beginning with XX. The 
number of "loop pairs" is recorded as Z, (M1365). 

To determine the "loops not touching each other, taken 3 ata 
time" (hereafter referred to as "triple products" or "loop triples”), 
it is now easier to use the loop numbers of the "loop pairs". These 


pairs are recorded as XX to Zoe 


The procedure is more readily explained by an example. Con- 


sider 5 loop pairs recorded as:* 


(1) (2) (3) (4) (5) 


XX AA AX, AX. XA, XX, XX XX, XX. XX o 


C B Bt+l D D+1 


The variables A, B, C, D are defined initially as the XX subscripts 
as indicated. In order to have any "triple products" one number of 
each of two pairs must be the same. By the manner in which the pairs 
were recorded, the first number of the pair is always less than the 
second and is always equal to or less than the first number of each 
Succeeding pair. Therefore, in order for loop number XX, to be 

one of a loop triple, XX, must equal XX. If not, then all four vari- 
ables are HCE SSIES by two and the 'check is made agin. If XXp = XXp, 
in order to have a.loop triple, loops XXo and XXpi1 must not touch 
each other, and hence will be recorded somewhere farther on as a 
pair. Check for this by using XXp- For this to occur XXp = XXa, if 
not then D is increased by two and checked again (until D = 2Z,). If 
XXpH = XXa, then to have a loop triple XX pi must equal XXpi1: If 
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this is found to be true, the loop numbers of the three non-touching 
loops are XX ps AX a and AAs These three numbers are record- 
ed as subscripted XX variables beginning with XX5> . The number 
of triple products is recorded as Zo (M1570), 

The "loops not touching each other, taken four at a time" 


' or "products of four'') are found by 


(hereafter referred to as ''4's' 
an extension of the above process, 


Consider the same 95 loops pairs as above, and two triple pro- 


ducts: 
(1) (2) (3) (4) (5) 
XX XX, XX XX XX, XX, XX XXy XX XXo 
A C B D 
(1) (2) 
RAig 241, S%r0 ®Xi3 SA, S55 
i 
" " 1 = = 
For a product of four’ to exist XX XX, XX and XX a XXpi1 
XXpat do not touch each other, If the latter statement is true, then 
XX as XXpip XXH a1 have been recorded as a triple product and the 


program follows a procedure similar to that described for the triple 
products. If any "4's" are found, the loop numbers are XXy, XX, 
XXpay AX ap and these are recorded as XX variables immediately 
following the loop triples. The first loop number of the ''4!s' 


be XX 


" would 


Pe + 3Z.," The number of "4's" is recorded as Z, (M1795). 
3.7. Loops and Products of Loops not Touching Paths (M1805-2500) 
The same reasoning regarding the extent of loop products rou- 
tines applies equally to this subroutine and the physical limitation is 
placed at having no more than products of three loops not touching one 
another and not touching a path. This is not considered a limiting fac- 
tor since the complexity of the system with greater than this number 
of loop products not touching paths would be such that other limitations 


previously mentioned would prevail. 
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The method used here is similar to that of the L.. Subroutine. 
The procedure is to find all loops not touching path 1, then all loop 
pairs and loop triples not touching path 1. After these are stored, 
the process Starts again with path 2 and so on until completion. 


Consider the following as having been previously recorded: 


2 direct paths; 4 loops 
(1) (2) (1) (2) (3) (4) 
sia laa 7 Be OE, oo ae 
Se. . 
YY ny YY p4179 Z2Z,7=0 ZZ .79 227,70 22,70 
3 pairs 1 triple no ''4''s"' 
loops loops 
1 2 1,2,4 
Tee 
2,4 


The pairs and triples are stored, as noted before, as : 


XX AX, XXo AB AA, XA, AX AX, AA 


1 2 1 + 2 4 1 2 4 


and the following variables have been previously defined: 


Yo number of paths = 2 
Zo number of loops = 4 
Z, number of pairs = 3 
Zo number of triples = 1 
Z, number of Nats" = 0 
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G, first ZZ number of loop 1 = 1 

Go first ZZ number of loop 2 = Gal 
Ge first ZZ number of loop 3 = Dtl 
Gy first ZZ number of loop 4 © Gl 


The procedure is to check the first number of path 1 against all numbers 
of loop 1, then the second nu mber of path 1, etc., until YY aay (0) is 
found, If a number is found in the loop equal to the check number of 
the path, the loop touches the path (M1855), the loop number is in- 
creased by 1 and the check number of the path reset to YY,- If all 
mumbers of the path are checked against a loop and no common point 
found, the loop does not touch the path and the loop number is record- 
ed as the next subscripted XX, in this example the first loop found not 
to touch path 1 is recorded as XX o (M1900). After path 1 has been 
checked against all loops (M1835), the number of loops not touching 
path 1 is recorded as i (M1925). In the example, say that loop 1, 3 
and 4 are found not to touch path 1. These are recorded as XXQ=1, 
XX. 8, XX, 174, Y,=3. 


To find the pairs of loops not touching path 1 and not touching 


10 


each other, it is necessary only to check the loop numbers of single 
loops not touching path 1, against the loop pairs, recorded as XX to 


ay This procedure is indential to that used in the loop products 


routine: 
(1) (2) (1) (4) (2) (4) (i), (2) (4) 
XX XX, XX, XX, XX, XX, XX, XX, XX. 
(1) (3) (4) 
XXq XX, XX XX,» (not yet recorded) 
Q Cc B E 


The variables A, B, E, Q are defined initially as the XX subscripts 
indicated and the recorded values of the XX variables are as indicated 


in parenthesis. The procedure is to check 7 against XX py if not 
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equal, increase A by 2 and check again. In the example eo) = XX. 


The variable C then is defined as Q+l1 (=10) and XX_, is checked against 


E 
XX qui: These are not equal hence loops (1) and (2) are not a "pair 
not touching path 1". By the previously described numerical order of 


recording the XX variables, it is seen that if XX is greater than 
ca CT 
next pair is checked in the same manner. Again ey AX y (A now 

equals 2) and XX ee (M2025) and found unequal. 
Now, however, XX is less than XX ay so the possibility remains 

that XXo and XX. are a pair not touching this path. C is then increas- 
ed by one (M2045) (to 11) and XX 
This time the two are found to be equal and a) and XX a are record- 


ed as XX... and XX Now A is increased by 2 (to 4) and XXo is 


then further checks against XX, are necessary (M2035). The 


is checked against XX 


is again checked against AX) ||. 


2 jee he 
checked against XX, (M1965), found to be unequal and furthermore, 


AX, is greater than er so*no further test is necessary using BO 
The procedure continues and no more pairs are found. The number of 
pairs of loops not touching path 1 as recorded as Y, (M2100). Since 
only one pair has been recorded there can be no "triples not touching 
path 1'"'". If more than one pair had been found, these would have been 
checked for triples in a similar manner to that used for loops and the 
number of triples is recorded as Y3 (M2335), If none are found, Y. 
is recorded as zero. 

The next step is to determine the end of the present path in the 
YY list (M2385, 2390), check to see if there are any more paths (M2400) 
and if so, repeat the process for the next path. 

The entire XX list is summarized in the next paragraph. In ad- 
dition to the XX variables, the following have been defined by this sub- 
routine: 


Number of loops not touching path 1 
Number of pairs of loops not touching path 1 


Y 
Xi 
XY Number of loop triples not touching path 1 
Y 
x 


5 { Same as Yy> Yo: Y2 applied to path 2 


3G = 





<a a eh Fe a eet 


x 
3Y 972 


*3Y 71 {Same as above applied to the last path 


at 
3Y 5 
3.8. Performance Function Subroutine (M2500-3560) 
Before discussing this subroutine, it is best to summarize the 


XX variables recorded in the routines of paragraphs 3.6., and 3.7: 


¥, 


XX, -" ctr cc cc ee BoopeRairs:= - 1-2 a 
. —_aa 
ees i Loop Triples - - - XX 7t3Zo-1 
XxX ------- Loop "4's" - - - - XX 
2Z,+3 Zo 2Z,+32,+42, l 
Define E = 2Z,+3Z25t42,, 
XXy, = —— > Legpeenot touching P, ~~ AAR+Y -1 
oi rane - - Loop pairs not touching P, - = MARTY +2Y5-1 
=—ry 2%, not touching i --- AAHTY +2Yy43Yq71 
Define F=E + ¥,t2Yot3Yo 
XX. --- = Loeepe net touching Po 
etc. for all paths. 
The last number recorded is XXoq: where G is defined as: 
Gs 2Z,t38LZot4ZotY t2Yot3YntYy +=-+- rE 


As stated in Chapter 2., this routine is written to provide a per- 


formance function of the form 


17273 479 


1 + GG, G, + GgG, + - - - 


G,G,G, * G,G. t= = - 
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The results of the routine are the subscripts of the G's with appro- 
priate algebraic signs. This routine is by-passed if not desired 
(M2501). The by-pass is effected as stated in paragraph 3.2., by 
entering a zero for the variable SS on the first data card. 
The Subroutine is quite naturally divided into two parts, namely, 
a Numerator and Denominator. Recall that the numerator is: 
=zP, [1+2L, 


ue 2 OE a 2 2 ae 


k k k 


In order to punch this function, it is necessary only to know 
where the particular values are stored in the computer, The method 
is very Straight forward and will be described in detail only where 
particular sections of the program might require clarification. 

The YY and ZZ lists, of path and loop components respectively, 
also contain the Summers, which are not desired in the equation. 
Hence, any number in these lists greater than the number of compon- 
ents is not punched (M2600), 

For convenience, the first YY number of each path is recorded 
(M2515 and M2635) in the Same manner as the first ZZ number of 
each loop was recorded earlier. All components (excluding summers) 
of path 1 are punched successively (M2600, 2605) with no sign of posi- 
tive and with a minus sign of the term of negative (M2585, 2595). The 
second term of the numerator (P,. yD Lj) = P, Li, Sg P, Lo wm . iG 


a Series of Y, terms, the first of which contains all the components 


of path 1 and ai the components of the first loop not touching path 1, 
This loop number is found in the XX list, the first ZZ number of the 
loop is known so the components of the loop are located and the entire 
term is punched, 

Throughout the routine, there are many occassions to use a par- 
ticular group of equations or commands, equations locating the number 
in the XX list, punch plus and punch minus commands, etc. To obviate 
the necessity of writing these Same equations many times, two "routine 
indicators" P and R, are uSed to transfer control back to the desired 


program locations, 


457 « 





Continuing the numerator, the next series of terms is 


2 ale 
2 


rk 
and the procedure is similar to that used before. After all terms 
have been punched, the same process is carried out for path 2 etc. 
until the entire numerator has been punched. To fully understand the 
mechanics of the program requires a step by step analysis of a sam- 
ple system. This is done for a sample system of Chapter 4. 

After all numerator terms have been punched, the program 
shifts to the denominator (M3083), Hence, the routine is simpler 


since the denominator of the performance function equation is 
rt2ok + 20, + 20 besa. 
Si So © o RE ni Mg 


The procedure, is the same as for the loops of the numerator and the 
components of each loop, loop pair, etc. are punched successively 


with appropriate algebraic signs. 


Sig wl Frequency Response Calculation (M3570-M5285) 

Having set up all terms of the performance function, there re- 
mains only the task of devising a scheme for entering the transfer 
functions of the system components, substitute 5 = jw, assign a fre- 
quency value and combine the results into the terms of the perform- 
ance function. As described in Chapter 2, the component transfer 
functions are a ratio, the numerator and denominator of which are 
comprised of products of terms of the form AS? + BS = C,. of wien 
any of the coefficients may be zero. It is desired to enter on data cards 
only the coefficients. Since each transfer function may contain any 
number of terms, provisions must be made for recognizing the number 
of terms for each component. This is accomplished by means of an 
array of the subscripted variable G. For ease of entering the infor- 
mation a complete form is included with the sample problem of Chapter 
4. 
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To understand the form in which the data is entered and the re- 
quired information obtained, a simple example is shown in Fig. 3.2. 
As stated previously, in order for the program to be properly 
executed, a scheme must be established whereby the terms of the 
array are recognized, The scheme chosen to accomplish this is to 
enter the G subscript of the first number of the last term of each 
numerator and denominator. These numbers will be referred to as 
"numerator stops" and ''denominator stops" respectively. The num- 
bers are shown circled in Fig. 3.2. and are entered in the following 
order: 
(1) Numerator stops listed consecutively, 6 to a card 
(2) Denominator stops listed consecutively, 6 per card, 
The last card of each may not necessarily be full. In the example of 
Fig. 3.2., only one card is needed for each. The first card would 
contain: 0, 15, 24 and the second card: 6, 21, 27. These numbers 
are entered into the G-array immediately following the last number 
of the transfer function array, with a blank space between the two 
sets of numbers. This blank space is given the value of -3, which 


is necessary for manipulations later on. 


Example of the Transfer Function Array 


Consider a system of 3 components with transfer functions: 


" eel 
S (58S + 1) 
2) 62 (,.2S + 1) (3S # 1) 
(ot ce + 7S 4 1) 
9 <2 
| 


The coefficients are entered in the following order: 
1. All numerator terms of component 1 

All denominator terms of component 1 

All numerator terms of component 2 

All denominator terms of component 2 


All numerator terms of component 3 


nm OFF FP W BD 


All denominator terms of component 3 
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The numbers are entered consecutively, 6 per data card. These be- 
come a G-array as (numbers in parenthesis and circled are Sub- 


scripts of the variable G): 


Conerant Coeff ers octet” 
Numerator 1 1 (|) oe (2) 0 
Denominator 1 (3) 0 (4) 1 (5) 0 
1 (Teme) (8) 0 
Numerator 2 (9) 62 (10) 0 (iy 8 
{ (12) 1 (13) 652 (14) 0 
(15) 1 (156) es (17) 0 
Denominator 2 (18) 1 (19) 4 (20) 0 


Fig. 3.2. G-Array of Coefficients of Transfer Functions 


One other card is required to be entered and it contains the fre- 


quency information as: (M3570) 


Initial frequency (radians /Second) - - - - - VV 
Frequency increment desired - - - - - - - DW 
Final frequency (radians/second) - - - - - FW 


Since a MAC READ order for reading in arrays of numbers must con- 
tain some information as to the number of cards to be read, it 1s neces- 
Sary that the last G subscript of the transfer function array be entered, 
This number (shown circled in Fig. 3.2. as 29) is entered on the fre- 


quency data card and recognized by the computer as RR. (M3570). 
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Having read in the necessary data, the program picks up the 
transfer function of the first component. Letting S = j w and sub- 
stituting the initial frequency (VV), a magnitude and angle is com- 
puted for the component. Then the Same process is carred through 
for each component. Referring to Appendix B, card M3640, the pro- 
cedure is as follows: 

1. Define the variable H initially as the subscript of the 
G variable left blank between the numerator and de- 
nominator stops. | 

2. Establish the variables K and L as thelimits in the 
G-array of all terms (numerator and denominator) 
of component 1. 

Space available as subsScripted XX is used as temporary storage 
for the computation of the magnitude and angle of the components, The 
magnitude and angle of each term of component 1 is computed by M3690 
and M3730. The command DO TO is used to accomplish this and means 
to carry out the operation specified for various values of the variable 
K from its present value in steps of (3) to a final value L. The opera- 
tion is merely to determine the square root of the sum of the squares 
of the complex number representing each term of component number 1. 
The angle is found by M3730 as the ARCSIN of the imaginary part of 
the complex number divided by the magnitude. The ARCSIN function 
will give results betweent!!/ 2, hence, it is necessary to make the check 
of M3760 which states that if the real part of the complex number is 
negative, Subtract the angle from II, If the real part is positive, the 
operation is ignored. 


To see the manipulations more clearly, refer to Fig. 3.2. 


R = 29 
H = 33 

K = Gag + 3 = 0 
L = Gg, = 6 


This is seen to indicate 3 terms of component 1, and by M3690 and 
M3730, a magnitude and angle is computed for each of them. To ob- 


tain the magnitude and angle of the entire component, the numerator and 
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denominator terms are Separated, each manipulated as necessary, 
then the magnitudes divided and the angles subtracted to give the 
results. This is accomplished by M3780 to M3925. Referring again 
to Fig. 3.2., and starting with M3785: 

K = 0 


M=G 0 


Zi) + 1 
This means that there is only one term in the numerator and the 
numerator magnitude then is the magnitude already computed for that 
term. The numerator angle is found by adding the angles of all 
numerator terms (M3835, M3840) but in this case, there is only one. 

A similar procedure is carried through for the denominator 
terms: 

(M3855) K = Gogt3 = 3 established S as the G subscript 

of the first denominator term. The two denominator terms 
are multiplied by M3865 and 3870, in two steps, K = 3, 6. 
M3850 Cel 
M3870 C = (1) (XX 44) = XX i3 
M3870 C= (XX (XX ) 


E+3) E+6 
In a similar manner, the 2 denominator angles are added. (M3895, 
M3900) 

It is desired now to compute a magnitude and angle for the com- 
ponent as a subscripted variable G in the interests of space economy, 
since later subroutines use the XX variable. Therefore, F is com- 
puted as the next available G number, and the component magnitude 
is computed by dividing numerator by denominator and identified as 
Gp. 


numbers followed by all component phase angles. Hence, with Q 


It is desired to have all component magnitudes in successive G 


{ 
components the phase angle of component one is computed by subtract- 
ing denominator angle from numerator angle and identifying the result 
as CEiQ: 

The procedure now is to return to the beginning of this subroutine 
and carry out the above process on the remaining components in suc- 


cession, Since the component magnitude and phase angles will be used 
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later, it is convenient to define the variable J (M3970) as the last G 
number preceding the first component magnitude. Then to find the 
magnitude of component N, the location is known to be Gyan: 

It would be satisfactory now to Substitute these components 
magnitudes and phase angles into the performance function as punch- 
ed in paragraph 3.8., perform the required algebra and obtain the 
response of the system for this frequency. However, the perform- 
ance function is not retained as Such and calculation of the perform- 
ance function is deSired to be optional. Since the XX variables are 
stored as loop numbers, each loop number appearing many times, a 
more convenient and faster method is to compute the magnitude and 
phase angle of each loop. Then the loop response is used in comput- 
ing the response of the loop product terms. All of theSe are terms 
in the performance function and will be referred to herein as "term 
responses". 

The first step is to devise a storage scheme for these term 
responses so they may be easily combined. The denominator of the 
performance function is computed first and the storage scheme used 
is (M4005 - 4025): First available XX variable: XX 


E 
ee ~~ Seer! eke Se — 
Loop magnitudes Loop phase angles 
—<—1° °° EAC si ©i<s Nai 
Real part of all Imaginary part of 
term responses all term responses. 


The denominator is a Summation of term responses: 


(1 + a} ae) 


After computing these term responses and converting to rectangular 
form, all that is necessary is to add the real parts (XX 44 to XX) 


and the imaginary parts (XX OQ, toXXa7 4 Z,+Z 


il 0 + Z3), to convert 


2 
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the result to polar form and retain it until after the numerator res- 
ponse is computed. 

In computing loop magnitudes recall that the ZZ list of loop com- 
pOnents (paragraph 3.4.) has a sign determined by the number of minus 


signs encountered. For example, a single loop is recorded as: 


ZL, = =3 
aL, = 1 
ZZ. eZ 
ZZ, = 0 





As explained in Chapter 2 and from elementary feedback theory 
the sign of this loop is (+). Hence, an initial multiplier of -1 is used 
(M4045) and an initial phase angle of 0 (M4050). Component magnitudes 
and phase angles have been stored as Grin and Gy+A+Q respectively. 
The loop magnitude and phase angle is computed (M4055 - 4115) using 
only components (summers are excluded), with the appropriate sign, 
and recorded as described in the storage scheme above. It is required 
that the loop response be retained in both polar and rectangular form 
(M4130 - 4160). 

The response of the terms representing loop pairs is obtained 
by (M4210 - 4290) using the loop numbers of the pairs from XX storage 
beginning with XXp and multiplying magnitudes and adding phase angles 
of the loops. Since these loop patr responses are never multiplied in 
the performance function equation, they are converted to rectangular 
form (M4250 - 4260) and retained only in this form. 

By exactly the same method, the term responses for the loop 
"triples" and ''4's" are computed and retained in rectangular form 
(M4310 - 4500). The denominator response is then obtained by adding 
all real parts of the term responses, Stored consecutively as XXy 41 
to XX as and all the imaginary parts, XX ay to XX aims: Tisa 
variable used to increase the XX variable by one each time a new term 


was added, hence it can be used as a check to determine the end of the 
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list to be added. There are T-1 numbers in each list. The result of 
the additions, a complex number, iS converted to polar form (M4580- 
4630) and retained as denominator magnitude, (DM) and denominator 
phase angle (DA). 

Note, that the conversion from rectangular to polar form is a 
routine that is to be uSed more than once, So the routing identifier (X) 
is used (M4570). As noted previously, the ARCSIN function given only 
angles between + Lj so if the summation of the real part is negative 
the correct angle must be computed (M4605), 

Obtaining the numerator response is only slightly more com- 
plicated since it is a Summation of terms multiplied by a path function 
and then the products are summed: 


ae) | eee — La @ .-oeel 


k 
The procedure is to determine the path response (P,) for path 1, in 
polar form, then determine the term response of each of the terms 
indicated in the equation that are associated with this path, perform 
the addition indicated in the brackets and convert to polar form, This 
result is then multiplied by the path response, the final result convert- 
ed to rectangular form and stored. 

The procedure is then repeated for the second path, and so on. 
After all paths have been completed the fihal Summation, 2 P. [ ate eee ] 
is performed and the result converted to polar form. 

The loop responses previously computed are needed here but the 
other denominator terms are not, hence, the Same storage may be used 
for the numerator terms. In the program of Appendix B, this is only 


partially true, in that it was chosen to use only the space after XXaiz ; 


(B=CtHt Zo) 0 
oe ie oe te 
Loop magnitude Real part of all Imaginary part 
and phase angle denominator term of loop responses 

responses 

age) | a | os 
Real part of all numerator Imaginary part of all numerator 
term responses, term responses, 
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The magnitude and phase angle of the first path are computed 
(M4690- M4760) in the same manner as that of the loops except that 
the sign of a path function is minus if the path in the YY list has a 
minus sign. The path magnitude and phase angle are retained tempor- 
arily as V and Wrespectively. (M4720 - 4730). | 

The numerator terms associated with the path are located in XX 
Storage from the storage scheme Summarized in paragraph 3.8. The 
first group of these terms (= L,) is a Summation of loops for which 
the real and imaginary parts of their response have been previously 
computed, so need only be rewritten as the appropriate XX variable 
(M4770 - 4825). 

The second group of numerator terms (2 PL, L i? and the third 
group (21 L, i? are determined exactly as the > HL, 3 and 2TI L.. 
terms of the SD esibiniot (M4835 - 4915 and M4925 - 5015). 

Having determined the-response of all terms associated with 
this path, their real parts + 1 and imaginary parts are summed (M5030- 
59085). The routing identifier X = 1 is used and the program jumps to 
the routine for converting rectangular to polar form (M5095), 

Up to this point, the following have been obtained: v/ Ww [ P/ A 
for 1 + 2 L. 
this path is then computed (M5100 - 5120) and retained as an XX 


variable. If any more paths exist (M5130) the same routine is carried 


k ree Ne The response for all terms associated with 


out for each, 
After all paths have been completed, the real and imaginary parts 


are Summed (M5155 - 5195), The storage locations of these values are: 


Li xs 4 | ete. 


a { 
| XX AK B+3 *“Bi4 


AK B+1 **B+2 


ees 
first path second path 
Bach time a path was completed B was increased by 2 in preparation 


for working space of the next path. After all paths are completed, B 


is set back to its original value (M5156) to permit summing. 
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The program now jumps to the conversion of rectangular to polar 
form routine (M5210), the result of which gives the numerator magni- 
tude and numerator phase angle. The system magnitude ratifo (M) and 
phase angle (P) are found (M5215 - 5220), The final result is desired 
with the phase angle given as a positive angle between 0 and 360 de- 
grees, hence, the two equations (M5235 and 5255). The possibility 
exists that M might be zero or the DM might be zero (M would be in- 
finite). In either case, special punch orders are provided (M5212, 5275 
Macdo237, 5253). 

The following information is obtained from the punched result: 

M Magnitude Ratio 
Je Phase Angle in Degrees (0-360) 
AA Log M (MAC computes logs to the base e so it is 
converted to the base 10) 
20AA Magnitude, Ratio in d.b. 
VV Frequency (radians/sec) 


VV/ 2m Frequency (cps) 


3.10. Special. Features 

There are certain special Situations that are not strictly a part 
of any of the general Subroutines of the preceding paragraphs, but 
might require explanation. One such situation might arise when the 
system performance function is known in factored form, perhaps from 
a pole-zero configuration. In this case, the block diagram of the sys- 


tem-is: 


Component 1 may consist of any number of terms. There are no loops 
and only one path. The calculation of frequency response of the system 
is made treating the system as one component and the entire program 
up to and including the performance function routine is by-passed, This 
by-pass is effected by the presence of a zero on the first data card in 
the position designated as "number of summers" (M0035, 0036). The 
number of loops (Zo) is recorded as zero in the by-pass (M2565) and 


the variable E, used to establish storage locations in the XX list is 


ZAG, 





assigned a value of zero. 

The "input data array" is not entered. The first card must con- 
tain six numbers since the READ Order (M0035) says read six values. 
However, the only ones with any meaning in this case are the third, 
where a zero is entered for S, and the second where the number of 
components is entered as 1. The remainder of the data cards are 
entered (M3570, 3585, 39 95, 3605) exactly as explained in paragraph 
E. 9. 

A second special situation arises when, having computed a fre- 
quency response, it is desired to change the frequency range, fre- 
quency increment, or a transfer function of one of the existing com- 
ponents. In this case, it is desired that all of the loops, paths, and 
combinations leading to a performance function be retained, and only 
that part of the program pertaining to the calculation of frequency res- 
ponse be repeated. A simple procedure for handling this is to let the 
program run to completion, have new data cards available, and start 
again at M3570. The program is set up in such a manner that all 
that is required is to push the start switch and the program will Start 
at M3570. 

If in addition to the above, it is also desired to insert a new 
component or change the configuration of the existing circuit, it is 
still not necessary to read the entire load deck again. A “program 
start'’ may be entered directly into the computer in basic machine 
language. This entails setting up on an"index register", the address 
of the appropriate command in the computer. This is a minor task 
but to ensure proper operation, a detailed step-by-step procedure 


follows: 


Detailed Step-by-Step 'Re-run'' Procedure 


l. Prepare data cards as in original situation. 

2. Set 1698 in the ADDRESS SELECTION SWITCHES on the main 
console. These switches are plainly labeled. 

3. Set computer CONTROL SWITCH TO MANUAL, 


age 
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Press PROGRAM RESET key. 
Press TRANSFER key. 

set CONTROL SWITCH to RUN. 
Press PROGRAM START key. 
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CHAPTER 4 
SUMMARY OF PROCEDURES FOR USING THE PROGRAM 


4,1. General 

In Chapter 2, the general plan of the program was discussed 
with the intent of describing its application to any computer. Chapter 
3 covered the details of the program, of Appendix B, written in MAC 
language for the IBM 650 and‘its associated equipment. The mechanics 
of the program and the method used to obtain the frequency response 
have been adequately discussed. In this chapter, a Summary of the 
procedures required by the user will be outlined with application to a 
specific example. | 

It is noted that MAC uses plain language for the basic program 
but that this program is then assembled in "machine language" prior 
to operation. The "load deck'' referred to in Chapter 2 is a machine 
language version of the MAC program. Any other information entered 
in the form of data cards must be in a form recognized by the computer 
and independent of MAC. Data cards contain only 10-digit numbers, 
six or less to a card, and the numbers are in "floating point" form. 
Floating point is merely a method of placing the decimal point ina 
number always written with the first digit as non-zero. The first 8 
digits are the numerical value and the last two locate the decimal point. 
For a number with the decimal point at the extreme left (e. g. .62), 
the 9th and 10th digits are 50. For each digit that the decimal point 
is located to the right or left of the first non-zero digit of a number, 


1 is added or subtracted from 50, For example: 


62.0 is written in floating point as 62000000 52 


6, 20 as 62000000 51 
. 620 as 62000000 50 
. 062 as 62000000 49 
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4.2. Data Cards Required ‘ 

The data cards are placed in the card reader feed directly on 
top of the MAC or load deck and are read at the appropriate time by 
a READ order in the program. The number of data cards required 
will vary with the complexity of the physical system being analyzed. 
However, the type of data falls into six general groups. A new card 
must be punched at the beginning of each of theSe groups Since they 
are read in by different READ orders. The last card of each group 
‘is not necessarily filled. The six general groups and the data cards 
of each are asfollows: | 

(a). General Circuit Information (1 card) 

(1) Maximum number of inputs to any unit. .... K 


(2) Number of component boxes in the block 


Gl agen Of te Sy Stes oo oo QQ 
(3) Netmser ofeummers. ...... ... = 
—) Wumber of the imput unit ...... . 5. seen 
(5) Number of the output , unit ee 2 ee. DD 
(6) Is Performance Function desired? ... SS=0if No 
1 if Yes 


(b) Input Data Array (M cards) 
(where M = QQ +S +1, unless Kis greater than 4, in 
which case, the number of values on each card exceeds 
6 and additional cards must be punched) 
The "input data array'' described in paragraph 2.2. and 
3.3. is entered, one row per card (unless the row con- 
tains more than 6 numbers, in which case, an additional 


card is used for each multiple of 6). Each row of the 


array is Started ona new card. The last row of the array 


and hence the last card contains all zeros. 


(c) Frequency Data (1 card) 
(1) The lower limit of the frequency range over 
which the response iS deSired......... VV 
(2) The increment of frequency over the fre- 


quency range. . . S455 «3 smeteee a DW 
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(3) The upper limit of the frequency range... FW 
(4) The last subscript number of the G-array 
oftin@wierer PurctienS., 2. 6 i. se Ye RR 
Note that only four numbers on this card will be read. The last 
two Spaces will not and zeros may be entered in these spaces if de- 


sired: 


Ge The Transfer Function Array 

The number of cards here will vary with the complexity 
of the physical system. Numbers on these cards are 
punched in succession, six to a card, until all transfer 
function coefficients have been punched. As described 
in paragraph 3.10., the order of the coeffibients enter- 
ed for each term is: constant, coefficient of S, coef- 
ficient of a The order in which the terms are entered 
is: 

All numerator terms of component 1 

All denominator terms of component 1 

All numerator terms of component 2 

All denominator terms of component 2 

And so on for each component; the last card may 


contain 3 or 6 numbers, 


e, Numerator Stops 
A new card is started for this group, regardless of 
whether the last card of group (d) is filled. Numbers on 
these cards, as explained in paragraph 3.9., represents 
the Subscript of the G number in the transfer function 
array, containing the first coefficient (constant) of the 
last term of the numerator of each component transfer 
function, These numbers are punched consecutively, 6 
to a card, except that the last card may contain less than 


six, 
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(f) Denominator Stops 
A new card is started for the first numbers of this group. 
As explained in paragraph 3.9., numbers on these cards 
represent the Subscript of the G number in the transfer 
function array, containing the first coefficient (constant) 
of the last term of the denominator of each component 
transfer function. These numbers are punched consecu- 


tively, 6 to a card. 


4.3. Results of the Program 

As described in Chapter 3 and illustrated by the example in 
this chapter, the result of the performance function routine (if such 
is desired) is a series of groups of numbers representing the com- 
ponent making up a term in the performance function equation. These 
groups of numbers (terms of the equation) are preceded by the appro- 
priate algebraic Sign. 

Results of the frequency response calculation are given, for each 
frequency, aS a row of six numbers, under proper heading and re- 
presents: 

Magnitude ratio 

Phase angle (in degrees) 
Log of magnitude ratio 
Magnitude ratio in decibels 


Frequency (radians/sec) 


mo 82 0 DO ® 


Frequency (cycles/sec) 


4.4, Example 

Consider the system represented by the block diagram of Fig. 4.1. 
with transfer functions as indicated. This is a hypothetical system 
chosen of sufficient complexity to adequately test all parts of the pro- 
gram and to illustrate the advantage of digital computation over graphi- 
cal techniques. 


The first step is to list the information required by the single card 
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Output Inputs 


1 8 0 0 
2 9 0 0 
3 10 0 0 
4 al 0 0 
9 12 0 0 
6 13 0 0 
¢ 8 0 0 
8 -6 0 0 
9 2 0 
10 2 =) 7 
11 3 -4 =) 
12 + 0 0) 
13 D -6 8 
O . 0 0 0) 


Fig. 4,2. Input Data Array 
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2G OOOOOD0G6E 


2G QOOOQ009E 


Ts OQOQO000T 


TS QOQOOOOOT 


cS 
6% 
OS 


TS 
TS 


TS 


TS 


OOOO00EE 
OO00000€ 
OO00000E 
OOQQ000S 


OOQ0000T 
H6V0000R 


OOO00000T 


OO000009 


cS 


cS 
TS 
TS 
TS 
TS 
TS 
cS 


TS 
cS 


‘ES 


TS 
TS 


OO0000Le 


O00000%72 
OOOOOOOT 
OOOOO0O00T 
OOO0000T 
OOOD0000T 
OOQODO0O00T 
OOO0000T 


OOO0O0O00T 
OOO000LY 


O0000008 


OO00000S— 


OOQ0000L 


TS OOO00008 


SaquvD VLVa 


cS OOOSSOT c 


cS OOQOOOOST 


TS OOQ00009 
TS OQO000009—- 
TS OO00000%— 


TS OOOQOO00E- 
TS 00000002 


TS OOOO0009 


cS 


Ts 


OG 
TS 
TS 


OG 


TS 
TS 
TS 
TS 
TS 
TS 
TS 
cS 
cS 
cS 
cS 
TS 
TS 
TS 


yey aloes 


O00000CT 


OOQ00009 


O0000002 
OOOOO0OO0T 
OOOOOOOT 


OO000004 


OOQ0000¢ 
OO00000% 


OO00000E 
O000000¢ 


OOOOOOOT 


@0000009- 


OO0000008 
OOIOPV0E I 
oo0000ddT 
OOOOOOTT 
OOOOOOOT 
O0000006 
O0000008 
OOQO000L 


cS 
TS 
cS 


TS 
TS 
TS 
TS 


TS 
TS 
cS 
OG 


cS 
cS 
cS 
cS 
TS 
TS 
TS 
TS 
TS 
TS 
TS 
TS 
TS 
TS 


OO0000SY4 
OOO00000E 
O00000¢4 


OOOOOODOT 
OOO0000T 


OO000000S 
OOOOOODOT 


OOOOOOOT 
OOOOQOOOT 
OOOOOOST 
OOOOOO0OT 


OO0000ET 
OO0000CL 
OOOOOOTT 
OOOOOOOT 
OO0000006 
00000008 
OOO00000L 
OO000009 
OO00000S 
OOO000004 
OOO00000E 
O0000002 
OOOOO0O00T 
Oo000000€ 
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of group (a). Values for the variables are seen to be: 


K =3 
QQ = 7 
Ss 6 
CC =8 
DD = 6 
sso" 


Succeeding steps ‘in preparing data for data cards are: 

1. Write the input data array as in Fig. 4.2. 

2. List the frequency data. For this example, a frequency 
range of .1 to 6.0 radians per Second was chosen witha 
frequency increment of .4 radians/sec. Note, that the 
last item on this card does not pertain to frequency data > 
but is a requiredentry and placed on this card for conven- 
ience. Its value will be determined later. 

3. Write the transfer function array as in Fig. 4.3. Indicate 
in some manner (circles are used in Fig. 4.3) the numera- 
tor and denominator stops and the last G subscript of the 
array. 

Data cards are now punched in floating point form. These are illus- 
trated in Fig. 4.4. 

The program, followed by the data cards, is then placed in the 
card reader and the computer set in operation. Since the first data 
card indicated a performance function was desired, the first result 
is this performance function. For this fairly complicated system, the 
performance function is quite lengthy so is appended as Appendix C. 

The frequency response results are listed and plotted as shown 
in Appendix D. Note, that a plotter is available and when Set to the 
appropriate scales will produce a plot directly from the punched cards. 


The plotter was not used for this example. 
4.5, Summary and Recommendations 


This report has illustrated the basic simplicity of the plan and 


the feasibility of programming such a plan for digital computation.,; 
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Some salient features are: 

a. The limit on size and complexity of the system to 
be analyzed is that caused by capacity of the com- 
puter, more Specifically, the capacity of the partic- 
ular method chosen. 

b. The time required for analysis is: 

1. Read load deck - - 5 minutes 

2. Computation of all the terms of the performance 
equation. This will vary with the complexity of 
the syStem. For the system used in example 
number 1, the time was approximately 8 minutes, 

3. Computation and punching the performance function. 
Again this time will vary with the system, For 
example number 1 the time was approximately 4 
minutes, 

4. Computation of frequency response for the system 
of example 1 - - Slightly under 2 minutes per fre- 
quency point. 

As a future project it is recommended that: 

a. The program be refined to reduce computation time. 

b. The refined program be altered to make use of aux- 
iliary equipment such as tape to lessen present phys- 
ical limitations. 

Ce With greater Storage capacity, the performance func- 
tion routine be extended to enable computation of the 
performance function as a ratio of polynomials of S. 

d. Having accomplished (c), an investigation be made 
into the feasibility of extending the program to include 
factoring the polynomials and obtaining a time res- 


ponse. 
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APPENDIX A 


Program Summary No, 1, L,, Subroutine (M0075 - M0585) 


The purpose of this Summary is to explain the details of the 
program covering the method explained in paragraph 3.4. Using 
the G-array of paragraph 3.3., the first summer listed is located 
by knowing the number of components and the ''width" of the array. 
For the sample system of paragraph 3.3., this starting point is 
Go» the value of which is (6). As in paragraph 3.4., the G-array 
is searched by rows, starting with G, (encountering a zero indicates 
the end of the row) until the search number is found. It is located in 
G,. To find the associated output (first number of the row containing 
Gz) the equation: L = M INTEGER (Lim is used; M is the width of 
the array and L is the G subscript (5). A new value of L is computed 
then as (4) hence, the number wanted is G, = (2). Continuing as in 
paragraph 3.4., no more outputs of (6) are found. The (2) is used as 
a search number and is found to go to (5) and also negatively back to 
i) 

Before proceeding further, it is best to look at the XX-array. 
(The numbers in parenthesis indicate the subscript of xX locating 


the position in the array). 


XX-ARRAY 

(0) 

6 

(1) (2) 

6 2 

(3) (4) (5) 

6 2 9 

(6) (7) (8) 
=o 2 7 
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Note, that the (-) appears on the first number of the line rather than 
on the number that caused it. The command that finds the search 
number in the G-array (M0230 and M0240) checks for both positive 
and negative numbers and if negative sets a variable R = 1 (0255). 
Then when the line is recorded, each term is multiplied by R (M0350). 
After the first number of the line is computed, Ris set to+l. This is 
desirable sihce in later computations, it is the sign of a loop that is 
important , so the sign of the finst number only is checked. 

The method of this routine is adequately described in paragraph 
3.4. The following is concerned primarily with the mechanics of the 
program. As noted earlier, before recording a line in the XX-array 
certain checks must be made, the first of which is: 'Check of a com- 
plete loop (M0265)". This checks the next number to be recorded 
against the first number of the line and if they are the same, the num- 
bers of that line represent units comprising a complete loop and are 
placed in the ZZ list (M0500). The next check is to See if the point 
is already recorded in the appropriate line. The newly found number 
must be checked against all numbers of the line (M0275 to M0285). 

If it is a duplicate point, a new line is not recorded but the search is 
continued. If not a duplicate point, then it must be checked against 
all summers previously used as starting points; these are recorded 
by (M0485) and the check is made by (M0320). 

The mechanics of recording lines in the XX- array igetairly con 
plicated and is most easily understood after getting farther along in the 
array, although the same commands are used to record the first few 
lines. Consider the stage of operation in paragraph 3.4., at which is 


recorded: 


(9) (LO) (ay ae 


lina") -6 2 7 3 (XXq to XX 5) 
(13) (14) (15) (16) 

Lieut) ae 7 4 (XX, to XX,,) 
(17) (18) (19) “Ze, 

line (g) -6 2 5 8 (XX 7 to XX oo) 


so = 





line (h) 


(4) 


line (i) 


(5) line (j) 


line (k) 


(21) (22) (23) (24) (25) (XX. to XX...) 
ea 5 8 1 21 29 


(26) (27) (28) (29) (30) 


-§ 2 ‘i 4 8 (XX, gto XX oo) 
(31) (32) (33) (34) (35) (36) 
-~6 2 ii 4 8 1 (XX, to XX qe) 


The array is separated into "blocks" by the solid lines and the blocks 


are given numbers indicated by the number to the left of the array. 


To handle this recording procedure, the following variables 


are used: 
W 
Y 


V 


FR 
N 


Number of current block (5) 

The XX subscript of the first number of the 
last line of the previous block 

The XX subscript of the first number ofa 
line containing the present search number. 
Present Search number 

In general, the XX subscript of the first 


number of the last recorded line. 


From the array above, assume line (j) has just been recorded using 


FF = XXo 4 = (8) as a search number; then 
W = 5 
V = 26 
N = 31 
Y4 = 17 
Y5 = 26 


Returning to search, no more outputs of (8) are found and it is desired 


to pick up XXo5 = (1) as the new value of FF. 


M0400 


FF = XX 


V=V-W.=a21 


V+We 1 “Soeeenee 
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Return to search the G-array and find (-1) going to (7) and after 
making the 3 checks noted previously increase NtoN+W +1 = 37 
and record (M0350) line (k) as 


Tone, 5, 8, 1, 7 (XX to XX 49) 
No more outputs from (1) are found and a new search number must 
be picked up. From the array above, it is seen that it is not desired 
to go back another line but to start a new block. The reason for the 
Yi variable is now apparent. Change V, as before, to V= V- W = 16 
and check it against Yy,_, to see if block 4 has been finished (M0405). 


V=Ne= 37 7 
BR = AX. y= AXy0 = 7 
W=Wt+i1==6 
Ye = 37 
N=N- 1 (The reason for this is not apparent at 


present but will be explained shortly) 


The variables are now Set up to go back to search and complete anr 


other block. The process continues until the XX-array appears as 





ie Ling th) -6 2 8 1 (XX,,_95) 
line (i) -6 2 8 (XX 630) 
i aie Cy) es 2 oF 1 (XX, _9¢) 
line (k) +625 817 OP gee) 
- im (1) 6-2 "5 2ga le me (XX 49_49) 
line (m) 625 8 17 4 (XX,5_c¢) 


YI 
iy 


Note that nothing from 
line (j) is recorded here 
since (1) goes to (6) and 
forms a complete loop; 
also goes to (7) which is 
a duplicate point. 
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Having recorded line (m) and set up to start a new block, the values 


of the variables are: 


V=WN = 50 
RR = 4 
W=Wt1=T7 
Y7 = 90; Ye = 3% 


N=N-1 = 49 


Return to search with FF = 4 and find it goes only to (8), a duplicate 
point. 
V=V- W = 43 


fy We SS 


Return to search and find (3) goes only to (6), a complete loop, so 
nothing is recorded. y 

V=V- W = 36 

V is less than Y¢ hence it is desired to start a new block 
but the routine, for this summer, is finished since there is nothing 
recorded in block (7). Recall that if anything had been recorded in 
block (7), N would have ben changed toN+W+1=57. Now\when 
V is set equal to N = 49, the check against Yw indicates the end. 
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APPENDIX C 


The numbers on the preceding pages of this Appendix, represent 
the component numbers of each term of the performance function of 
the sample system of Chapter 4. To illustrate the complexity of the 
system more clearly, the performance function is written in the 
more conventional form below. Since there are so many terms the 


numerator and denominator are written separately. 


NUMERATOR: 
Gg - Gg Go + Gg Gg + Gp Gy + Gg Gy Ge - Gg Gy Gg ~ GeiG Grae 


-Ge Go Gy Ge + Ge Ge G, + Gg Ge Gy G, - Ge Go Ge Gy - Gg Go Gz G, G, 
+ G7 Gg Gy Ge Gg - G7 Gg Gy G5 Gg Go + G, Go Gg G, Gs Ge 
DENOMINATOR: 

1+ Gep+G, Gg Gy Ge Gg tG, Gy Gg Gy Gs Gg - Go + Gg + G, + Gy Ge + 
+ Gp - Gg Go + Gp Gg + Gp G, + Gg Gy Gs - G7 Gg Gy G5 Gg Gy - 

~ Gy Gg - Gy Gy - Go Gy Ge - Go Gg +t Gg Gy + Gg G, G, + Gg Gp + 

+ G4 Gg t+ Gy Go Gg - Gg Go Gg - Gg Go Gy - Gg Go Gy Ge + Ge Gg Gy + 
+ Gg Gg Gy Gs ~ Gy Gg Gy - Go Gg Gy Gg - Go Gg Ge - Go Gy Ge - 

- Gp Gy Ge Gp + Ge, Gy Gg + G, Gy Ge Gp - Ge Go Gg Gy ~- Gg Gy Gg G, Ge 
~ Gy Gg Gy Ge ~ Go Gg Gy Gs Gg 
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